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We study, both analytically and numerically, an ARCH-like, multiscale model of volatility, which 
assumes that the volatility is governed by the observed past price changes over different time scales. 
With a power-law distribution of time horizons, we obtain a model that captures most stylized facts 
of financial time series: Student-like distribution of returns with a power-law tail, long-memory 
of the volatility, slow convergence of the distribution of returns towards the Gaussian distribution, 
multifractality and anomalous volatility relaxation after shocks. At variance with recent multifractal 
. models that are strictly time reversal invariant, the model also reproduces the time asymmetry of 

Oh' financial time series: past large scale volatility influence future small scale volatility. In order to 

quantitatively reproduce all empirical observations, the parameters must be chosen such that the 
q ■ model is close to an instability, meaning that (a) the feedback effect is important and substantially 

ryj ' increases the volatility, and (b) that the model is intrinsically difficult to calibrate because of the 

^ | very long range nature of the correlations. By imposing consistency of the model predictions with 

^ . a large set of different empirical observations, a reasonable range of the parameters value can be 

determined. The model can easily be generalized to account for jumps, skewness and multiasset 
correlations. 
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INTRODUCTION 



The quest for a faithful mathematical model of price fluctuations has been taunting researchers for more than a 
century now, starting with Bachelier's random walk model in 1900 1]. Such an endeavour is important for a bevy of 
, reasons, both from the point of view of (a) fundamental economics (what is the cause of price variations and what 
information do they reveal?) and (b) of financial engineering, with option pricing, risk control and trading models as 
obvious applications. 

In an ideal world, "the" mathematical model of price changes should be simple enough to allow easy calculations 
and calibration, yet rich enough to embrace all known stylized facts that the recent access to huge amounts of data 
has helped establish. It is now widely accepted that price changes reveal (i) fat tails, well described by a power-law 
decay of the probability distribution for large returns Q, |& y|, (ii) long range memory in volatility fluctuations or 
volatility "clustering" , again described by a power-law decay (in time) of the autocorrelation of the volatility H, IE 
i-£h . and (iii) asymmetric causal correlations between past price changes and future volatilities, often referred to as the 
"leverage effect" [g| (for reviews, see e.g. [tj El El ^3)- We discuss below other, somewhat related, stylized facts 
^ . that have been reported in the recent literature, such as multifractal scaling, critical relaxation of the volatility after 
a shock (the financial analogue of the Omori law for earthquakes), etc. More recently, some statistical asymmetry of 
financial time series under time reversal was pointed out |1 3| - in other words, financial time series do distinguish 



'_ 

past from future. This might appear trivial but constitutes in fact, as we discuss below, a very strong constraint on 
the family of eligible models for financial time series - for example, Bachelier's random walk model is strictly time 
reversal symmetric. 

Scores of different models have been proposed to improve upon the simple Brownian motion model, which has 
neither fat tails nor volatility clustering. Levy processes allow one to sup erimpose jumps to the Brownian motion, 
and therefore generate fat tails, but has no volatility clustering [l(J, HH LL-J 03 GARCH models or simple stochastic 
volatility models such as the Heston model allow one to get both fat tails and some sort of volatility clustering, but 
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not the long memory observed in the data pH UtI ITM IT^|. Models that mix jumps and stochastic volatility have 
been investigated l2Cl. Multifractal stochastic volatility mod els, initiated by Mandelbrot, Fisher and Calvet |2l| and 
much studied since^l IlllllilSElSElSEllllSIS, seem to capture in a parsimonious way a 
large amount of empirical properties. However, most multifractal models are again strictly time reversal symmetric 
and lack an intuitive interpretation in terms of agent based trading models |35j . We in fact strongly believe that any 
serious model of price fluctuations should in fine be justified by reasonable behavioral rules and market microstructure 
effects (see [M H3 S S E3 for recent work in that direction). Quite recently, one of us (LB) has proposed, in 
the context of option pricing, a "statistical feedback" process where the local volatility is large when price moves are 
deemed rare, leading to a non-linear diffusion equation for the price This equation can be solved and leads to a 
Student-Tsallis distribution for price changes at all times 42] . In its original form, however, the model breaks time 
translation symmetry: there is a well defined starting date and starting price. Although this can be used to price 
options [ll], [4|| (in the spirit of the Hull- White model for interest rates |44(), the process has to be modified to be 
interpreted as a bona fide model of returns. Such an extension, and its modification to account for long-range memory, 
was proposed i n lip a nd recovers, following a different route, a multiscale GARCH model proposed by Zumbach and 
Lynch in 2003 [jjiEil (see [affilEa] f° r earlier work in that direction). Numerical simulations of this model suggest 
a very rich phenomenology, that seems to account for most stylized facts of financial time series. 

The aim of the present paper is to motivate this new model, discuss its relation with previous work, and investigate 
in full details its statistical properties, both analytically and numerically. We focus in particular on the probability 
distribution of returns which is the crucial ingredient for option pricing and risk control. Although not an exact result, 
we find that these distributions can be well fitted by a Student-Tsallis form, with a lag-dependent tail exponent. We 
reproduce in great details most empirical facts, including the anomalous relaxation of the volatility after a shock, and 
the past/future asymmetry of the time series. The model can be generalized to include jumps, the leverage effect, and 
multi-stock correlations. We then discuss the issue of calibration. Within strict econometric standards, calibration is 
extremely difficult due to the long-memory nature of both the empirical volatility process and the theoretical models 
that are constructed precisely to capture this long memory. We advocate the idea of 'soft' calibration, which in such 
cases should consist in reproducing semi-quantitatively as many observables as possible. These observables should be 
chosen to be robust to the details of the model specification, and test different "orthogonal" predictions of the model 
(these statements will be made clearer in the course of the paper and in Section IVT) |. Consequences for option pricing 
are briefly discussed, and will be the subject of another paper. 



II. SET UP AND MOTIVATION OF THE MODEL 



In the following, we will consider a discrete time model, with an elementary time scale equal to r, for example r =1 
minute. [A continuous time version of the model will be discussed below]. The price at time ti = ir will be noted p^. 
We will conform to the standard of dealing with the log-price Xi = hipi and define returns as r, = i^+i — Xj . |59| The 
random return is constructed as the product of a time dependent volatility cr, and a random variable ^ of zero mean 
and unit variance: 

r< = (it + Uiii\[r, (1) 

where \i is the average drift, which we will set to zero in the sequel, meaning that we measure all returns relative to 
the average drift. The noise £j can a priori have any probability distribution to account for high frequency kurtosis 
and jumps, but for simplicity we will mostly focus in this paper on the case of a Gaussian noise. However, as we 
discuss below, the introduction of jumps is needed to faithfully reproduce real price time series. 

The seminal insight of ARCH or GARCH models |4j| is that the volatility process reflects trading activity and is 
subordinated to past price changes. Intuitively, the level of activity becomes high when past price changes are, in 
some sense, anomalous. In the simplest ARCH model, this is expressed as: 



^-2 ~2 



r 
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(2) 



meaning that the volatility is equal to its 'base level' <7q plus a contribution coming from the last price change. In 
fact, we have written the feedback term in a way that expresses the comparison between the square of the last return 
and its expected value, equal to (TqT. If the last price change was small compared to usual, the volatility today is 
close to its normal value, whereas in the other limit, the last return is deemed anomalous and leads to a potentially 
large increase of today's activity. 

An argument motivating Eq. above is as follows. Suppose that some traders open positions (for example, long) 
at time ti—±, when the price is Pi-\. Such trades are often initiated with both a profit objective and a risk limit, 
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FIG. 1: Schematic shape of the distribution of stop loss/stop gain thresholds around the opening price of the trade. Plain line: 
symmetric distribution; dotted line: asymmetric distribution, giving rise to the leverage effect discussed in section IVI 



which would close the position at time U if the price has moved up too much (stop gain) or down too much (stop 
loss). It is very natural to hypothesize that to each opening trade are associated two thresholds, one above, one below 
Pi-i, that trigger a closing trade if exceeded. If many agents open both long and short trades at U-±, one can expect 
a quasi-continuous distribution of thresholds at pi_i(l + A), more or less symmetrically distributed around pi-i, 
triggering with equal probability sell back or buy back orders. The density P(A) is, in the simplest case, even (but 
see Section El for the inclusion of the leverage effect) and obviously vanishes at A = since nobody opens a trade to 
close it immediately (see Fig. 1). The width of P(A) is given, in order of magnitude, by (Tq^/t since this gives the 
natural scale beyond which an event might be deemed anomalous. Hence, a relative change of price 7"j_i will trigger 
on the order of: [6(J 

iVi(n-0« I ' ll P(A)dA (3) 
Jo 

stop trades. These trades of random sign lead, on the next day, to an increase of the volatility as: 

of - ol + GNi/r, (4) 

where G is the average square impact per trade, and a 2 is the volatility due to all other trades. Taking into account 
that P(A) extends over a range aoy/r, one finally obtains a general single time scale ARCH model: 



a 2 



1 + G(^) 



(5) 



where the function Q depends on the detailed shape of -P(A). Taking for simplicity, in accordance with the above 
discussion, 

P(A) - P1 ^ CXP( -^' (6) 

(where is a number setting the width of the distribution of thresholds, and Pi the total number of opened trades) 
finally leads to: 

G(u) = 2g{3 2 (1 - exp(- U 2 /2/? 2 )) , (7) 

where g — GPi/2f3 2 a 2 T is the ratio measuring the impact of all stop trades compared to that of all other trades. The 
simplest ARCH model Eq. © corresponds to the limit u <§C j3, that is, neglects saturation effects related to the fact 
that stop limits are not placed arbitrarily far from the entry point (i.e. f3 is finite). When this saturation is neglected, 
G(u) is simply given by gu 2 (but see below, Fig. 11). 
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Although the above feedback mechanism is most probably at play in financial markets, a strong limitation of the 
above model is to consider that all traders have the same time horizon, equal to t in the above formulation. However, 
it is well documented that the activity of financial markets is fueled by traders with different time horizons, from a 
few hours to a few months or even years (see e.g. ^3,E3)- Therefore, stop losses or profit objectives are not placed 
only around the last price Pi-i but around possibly all past prices Pi-i, £ = 1,2, .... Correspondingly, the width of 
the distribution of these thresholds is calibrated to the volatility of the price over the particular trading horizon, i.e. 
(Jq^/Jt. The generalization of Eq. (J5J to this situation therefore reads: 



i+X>( 



\Xi - Xi_ e \ 



(8) 



Expanding Qi for small arguments finally leads to the symmetric version of the model studied in the present paper 
(the inclusion of asymmetry will be discussed in Section 0) : [6l| 



91- 



(xi - Xi-i) 2 



with 



Xi+l 



a%£r 



(9) 



(10) 



The coupling constant gg is proportional to the number of trades Pg with horizon £. Because traders with a longer 
horizon have slower trading frequencies and under-react compared to short term traders, it is reasonable to imagine 
that gg is a decaying function of £. Both for simplicity and because it allows us to reproduce several stylized empirical 
facts, we will choose gg to be an inverse power: 



gi = g/e° 



(ii) 



but other choices are possible. For example, Zumbach and Lynch have presented evidence that gi has additional 
peaks on the day, week and month times scales. These authors have proposed a model very close in spirit to Eq. JjJJl, 
and discussed some of its properties. In fact, Eq. © is a special case in the family of quadratic ARCH models, where 
the volatility is expressed as a general quadratic form of past returns: 



rjr k 



(12) 



which contains ARCH, GARCH, etc. Our specification insists that only combinations of returns 'reconstructing' 
actual price changes over different time scales occur in the above sum, because they correspond to quantities directly 
observable to the crowd of traders, which, we argue, strongly influence the trading at time i. Our model corresponds 
to a particular choice for M above: 



M(i;j,k)= f< 

i— max(i— — k) 



(13) 



whereas most ARCH models correspond a certain regression on past instantaneous square returns, i.e., to k) = 

K(i — j)5jk, with a certain kernel function K , usually corresponding to an exponential moving average, K(£) — c/ . 

With a power-law specification for gg, and the choice of a Gaussian distribution for the noise term £ in the definition 
of returns, our model is fully determined by only four parameters: <7 sets the volatility scale, r sets the shortest 
time scale over which feedback effects are effective, g measures the strength of these feedback effects and a describes 
the relative importance of short term traders and long time traders in the feedback process. It may however be that 
the assumption of a Gaussian noise for £ is insufficient to account for the high frequency statistics of the returns. 
In particular, one expects that true 'jumps' related to unexpected news are not described in terms of a volatility 
feedback process. It is easy to extend the model in that direction and choose another distribution for £. In the 
following sections, we will present several analytical and numerical results of the Gaussian version of this model, 
and compare them to empirically known results. But before doing so, let us give the continuous time formulation 
of the same model, which can be convenient for some applications, such as option pricing. Introducing the standard 
Brownian noise dW t , one may write: 



dx t = o- t dWt, 



(14) 
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with: 

* t =ob +ff r y^dt (t _ t/ + T)1+a - (15) 
This model is well defined as soon as a > 1, which is the case we will focus on in the sequel. 



III. ANALYTICAL RESULTS 



A. Unconditional distribution of the volatility 

Although our model (Eq. ©) expresses the volatility as a deterministic function of the past prices and the only 
source of randomness comes from the noise £j in Eq. (JTJ, the volatility effectively appears as a random variable, and 
one can ask questions about its distribution, correlations, etc. The simplest question concerns the average value of 
the volatility, which also coincides, for a stationary process, with the long term volatility of the price. Averaging will 
always be denoted below with brackets (. . .) around the quantity which is averaged. For the average volatility, one 
has (assuming stationarity) : 



oo w 00 

k 2 > - ^ 2 +E^ 7 = -o 2 + E*]<o*>. ( 16 ) 

t=i ' ' i=i 

This equation has a well behaved solution only if: 

oo 

z 2 = ^gi<l, (17) 

i=i 

where the above equation defines z%, the subscript '2' refers to the fact that we study here the second moment of the 
volatility. When zi < 1, the square volatility is amplified by a factor 1/(1 — z%) compared to the initial value ctq. In 
the case > 1, on the other hand, the process becomes non stationary and the volatility grows without bound as 
time elapses. It is clear that the condition Z2 < 1 can only be met if the sum of gi converges, which imposes that 
the exponent a is larger than one. For a > 1, one finds Z2 — gC{a), which delimits a region in the plane g, a where 
the process is stationary. In the following, we will often assume that a is larger than unity but close to it (which 
is suggested by empirical data), and use in this limit a continuous approximation for discrete sums. In particular, 
C(a) w l/(a — 1). We will find below that empirical data on stocks favors values of Z2 ~ 0.85 — 0.9, meaning that the 
square volatility is increased by a factor ~ 6 — 10 compared to its initial value Ctq. Therefore feedback effects might 
be an important cause of the excess volatility in financial markets ED • 

In order to compute higher moments of the volatility, one needs in general to know the full temporal correlation of 
the volatility, that we will establish in the next paragraph. Simplified, approximate calculations can be performed in 
two extreme cases: (i) no temporal correlations (ii) full temporal correlations. This leads to an equation for (er 4 ) of 
the form (1 — Z4)(er 4 ) =rhs, where the right hand side is finite whenever a > 1, and can be computed if necessary 
(see below). The important discussion concerns the value of Z4. We will denote Mk = A4(0; —k, —k), which behaves, 
for large k, as k~ a /(a — 1). Using the results established below (see Eq. l|2"o|l1. one can obtain a lower bound 
and an upper bound z^ 7> on the value of Z4. If correlations arc neglected, one finds: 

00 

z 4 >z 4 , < =3g 2 Y,M 2 k . (18) 

k=l 

If on the other hand, if correlations are overestimated and taken to be constant in time, one finds an upper bound for 
Z4 that reads: 

(00 00 00 \ 

E Mk} 2 + 2 e M l + 4 E( fc - v* 1 * ■ ( 19 ) 
k=l k=l k=l I 

As long as Z4 < 1, the fourth moment of a is finite, but if Z4 reaches unity, it does diverge, leading to an infinite 
kurtosis for the returns. For a — 1.15, we find 2:4 < = 0.16 z\ and Z4 > = 1.44 z\. This shows that the kurtosis k is 
certainly finite for Z2 < 0.833; numerical simulations below suggest that k indeed remains finite beyond that value. 
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The above argument is easily generalized to higher even moments of a, leading to an equation (1 — zm){<J 2n ) =RHS 
with: 

oo 

z 2n , < = (2n-l)\\g n J2 M k' (20) 
fe=i 

and a more cumbersome expression for z 2n ,>- For large n, one finds Z2n,< ~ (2gn/e) n , showing that however small 
the value of g, sufficiently high moments of the volatility are divergent. Since ti — a^i, the even moments of the 
returns are given by: 

(r 2n ) = (2n- l)!!(cr 2 "); (21) 

therefore high moments of returns themselves diverge, suggesting that both the unconditional distribution of volatility 
and returns have a power-law tail (possibly multiplied by a slow function), with an exponent equal to the order of the 
last finite moment. We will confirm this prediction numerically in the following section. Remember however that the 
above discussion is only valid when the noise £ is Gaussian; if £ itself has a non zero kurtosis, then its contribution 
should be taken into account. 



B. Temporal correlations of the volatility 

A well known stylized fact is that the volatility is a 'long-memory' process, which means that the temporal correla- 
tions of the square volatility decay as an inverse power of the time lag, l~ v , with an exponent v less than unity. This 
property turns out to be extremely important because it is at the root of the very slow convergence of the distribution 
of aggregated returns towards the Gaussian. More precisely, the kurtosis of the return ar» — over scale £, itself 
decays as £~ v instead of £ , which is the case when the volatility process has a short memory. Since the empirical 
value of v is, for stocks, on the order of v = 0.2 — 0.3, the slowing down is substantial and essential to explain why 
long dated options still have a smile. 

We therefore turn to the calculation of the correlation function of the volatility, defined as: 

= I- (22) 

In the limit g <C 1, one can quite easily perform a perturbative analysis that neglects terms of order g , to get: 



T{£) = 2g 2 



,^ fc 1 -" x ^ j 2 

^ (£+j) 1+a + ^ k 1+a (£+j) 1+a 

0<k<j y •" 0<j<k y •" 



(23) 



An analysis of this result for I 3> 1 finally gives, for a > 1 but close enough to unity such that one can use continuous 
integrals instead of discrete sums: 

orl (a) 

leading to a kurtosis exponent v = 2a — 2. The volatility is a long memory process whenever v < 1, i.e. 1 < a < 3/2. 
Comparison with empirical data, done below, suggests that a is in the range 1.1 — 1.2. The exact equation for !F{£), 
not restricted to small g 2 , can also be written down, although it is more cumbersome. For this calculation, one should 
note that averages such as (of£?) c (where the subscript c denotes a connected average) are non trivial, since the 
volatility randomness comes entirely from past returns themselves. This contrasts with many stochastic volatility 
models where the volatility Oi and the noise £j are often chosen to be independent (unless one wants to model the 
leverage effect). In the present case, one finds, for j < i: 

(^ 2 )c = g 2 M (°' ~ k > - k ')(vi-k<7i-k>Zi-kZi-k>$)c = 2g 2 (a 2 )M i _ j , (25) 

fe,fc'>0 

Now, the full self-consistent equation for T reads: 

oo oo 

T{£) = g 2 [3T(0) + 2}Y / M k M k+e + 4g 2 ]T M k M k+l [l + T(k - k') + 2g 2 M k ^] 

k=l k>k' = l 

oo oo £ 

+ 2g 2 ]T M k M k , +e {T(k - k') + 2g 2 M k - k ,}+ g 2 Y, E M k M k ,[T{£ ~ k' + k) + 2g 2 M^ k , +k ] (26) 

k>k' = l k=lk' = l 
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Specializing to I = leads to the fourth moment of the volatility studied in the above section. The two assumptions 
made there to obtain a lower and an upper bound correspond to Til) = ^(O)^^ and T(£) = J~(0), respectively. 
For large £, an asymptotic estimate of the various terms leads to the same decay as that predicted by the above 
perturbative calculation, i.e., J-(€) ~ J r 0O £^ l/ ) with v = 2a — 2, and a prefactor T^ increased by a factor 1/(1 — 
However, sub-dominant terms also appear, proportional to l~ 2v , £~ a , etc. The finite i behaviour oiT(€) would require 
to solve the above equation numerically. 

From the knowledge of Til) one can obtain the I dependence of the kurtosis of the returns, following Again, 
one should take care of the terms involving (of £|) c , which, as we discuss below, lead to a new, perhaps unexpected 
effect. One finds, for the kurtosis of the returns on lag t. 



k(€) 



k(1) 



■fl£a- 



-)[Hj) 



2g 2 M j ] 



For large lags £ ^> 1, one finds, using Eq. (|24|) . and for a close to 1: 



k(£) 



(3-2a)(2-a) 



(27) 



(28) 



Therefore, one expects the returns to converge to Gaussian, but only on a very long time scale. Any measure of the 
distance from a Gaussian - such as the mean absolute moment studied below - will tend to zero very slowly, as i~ v , 
see Figs 4-a, 4-b. If one now studies Eq. I|27|l for small values of t , say t = 2, one finds: 



k(2) - k(1) = 3[.F(1) - T(0) + 2g 2 M 1 



(29) 



In many models, the last term is absent, and since ^(l) < ^"(0), one usually finds that the kurtosis of aggregated 
returns is less than the kurtosis of elementary returns. However, the third term in the above expression suggests 
that one can observe, in some cases, a kurtosis that first increases with lag before decaying to zero. We will see that 
this is indeed the case in the numerical simulations of our model, although this effect is, again, very sensitive to the 
assumption that is a purely Gaussian noise. 



C. Conclusion 



The summary of this technical section is that the two major stylized facts (fat tails and volatility long-memory) 
are present in our model. We have indeed shown that the distribution of returns and of the volatility have power- 
law like tails, since high moments of these distributions diverge. We have also shown that the temporal correlation 
of the volatility is decaying as a slow power law. The following sections will be to establish these properties more 
quantitatively using numerical simulations, and to show that many more stylized facts can be reproduced by the 
model. Finally, we will turn to the question of calibration and discuss how the model parameters can be chosen to fit 
empirical data. 



IV. NUMERICAL RESULTS 



We have established above that the volatility-volatility correlation function, and the kurtosis, decay at long times 
as t~ v with v — 2(a — 1). A large amount of empirical work on financial time series suggest that v is in the range 
0.2 — 0.4 for many different assets. For example, averaging over the 500 largest stocks of the NYSE leads to v ss 0.25, 
while v w 0.3 for the S&P 500 Index [HI- We therefore choose to fix a = 1.15 (corresponding to v = 0.30) in most 
of the numerical experiments that we have conducted. Other values of a are briefly discussed, in particular in the 
context of the model calibration. The choice a = 1.15, although guided by empirical data, immediately leads to a 
numerical problem due to its proximity with the critical value a = 1 which separates a (theoretically) stationary 
regime for a > 1 from a non stationary regime for a < 1. The convergence of (say) the average volatility to its 
asymptotic value is expected to occur at speed T 1-Q , where T is the total length of the time series. For a — 1.15, 
this is extremely slow: even for T = 10 6 r, one expects corrections of order 10% to the theoretical asymptotic results. 
For this reason, and also to speed up the numerical calculation of the sum that determines the volatility (Eq. 0), we 
have truncated the power-law memory kernel gg beyond I = 5. 10 4 . The total length of our simulations is usually 10 6 
steps, but we discard the first 15. 10 4 points of the series before we start measuring any observable. Although this is, 
again, insufficient to obtain very precise results for such low values of a, we believe that these numerical experiments 
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FIG. 2: A typical time series (of length 2 10 5 ) of the volatility a 2 , for Z2 = 0.85 and a = 1.15. We have in fact shown a 300 r 
moving average of <xf, aimed at representing the 'daily' volatility within our model. 



are sufficient to obtain a good estimate of a host of different interesting observables, in any case comparable in quality 
to the corresponding estimates on real price time series. As will be clear below, we estimate that one day corresponds 
in our model to £ ~ 300; therefore 10 6 time steps corresponds to 3000 trading days, or twelve years of data. In the 
following, the base volatility ao is set to ao = 1, any other value would only change the following results by a trivial 
multiplicative factor on the returns. We will vary the coupling constant g, which we will in fact express in terms of 
z 2 = X^ff^> smce we know from the above discussion in section ITTT1 that it is really z 2 that measures the strength of 
feedback effects on the volatility. In the limit z 2 — > 1, we know from section ITTll that the volatility will blow up and 
the process becomes non-stationary for all values of a. Therefore, studying numerically values of z 2 too close to unity 
will also be difficult (the convergence is now as slow as [(1 — z 2 )T] 1 ~ a l), but, ironically, corresponds to the empirical 
situation. In the following, we restrict our simulations to the range z 2 € [0.60,0.85] - smaller values of z 2 lead to 
a process which is only weakly non Gaussian, whereas larger values of z 2 give rise to a numerically very unstable 
process, even though in theory the process should still be stationary on extremely long time scales. We will see below 
that values of z 2 as high as 0.9 might be needed to fit the data, but we have not attempted to simulate the model for 
such a large value. 

Although the issue of calibration will be more deeply discussed in section IVTl we will compare in this section our 
numerical results to empirical data, averaged over a set of 252 US stocks, chosen among the most liquid ones, during 
a four year time period: 2000-2003. 



A. Volatility distribution and volatility correlations 

1. Volatility distribution 

We first focus on the properties of the 'true' volatility <7j, which we can of course measure numerically but is 
unobservable directly in practice: only proxies of the volatility, obtained by averaging over several time steps, can be 
studied. A typical time series of a 2 is shown in Fig. 2, and reveals apparent shocks and volatility clustering familiar 
in financial time series. We show in Fig. 3 the histogram of u = lner for different values of z 2 . Obviously, since 
a > o"o = 1; the probability distribution function (pdf) of u is zero when u < 0. We have found that the pdf P(u) of 
u can be very accurately fitted by the following form (see Fig. 3): 



P{u) = Zexp 



u 



e(«), 



(30) 



where &(u > 0) = 1 and Q(u < 0) = 0. We have no detailed justification for this specific functional form for 
u — > 0. On the other hand, it is easy to show that the exponential tail for large positive u translates into a power-law 



distribution for a itself, decaying as a 



which is indeed expected from our theoretical analysis. Correspondingly, 
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FIG. 3: Histogram of u = lncr for 22 = 0.85 and a = 1.15, and two fits, using Eq. iJUJ - "Fit f3" - and Eq. (JSU - "Fit 
Student". We also show the slope fj, = —3.5 for comparison with the tail of the distribution of u. The pluses correspond to the 
histogram of the average volatility over 100 time steps, close to what one would determine empirically from price time series. 
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0.69 


3.1 
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0.54 


3.26 


0.85 


3.7 


0.80 


5.0 
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0.51 


5.67 


1.03 


4.2 


0.85 


6.66 


6.05 


3.95 


4.03 


0.50 


6.83 


1.24 


5.5 



TABLE I: Value of different observables and fit parameters for different values of 22, at fixed memory kernel (a = 1.15, or 
v — 0.3). The values of A 2 must be multiplied by 10 -2 . Note the 10% discrepancy between the theoretical and empirical value 
of (a 2 ) when 22 reaches 0.85. The value of ^ suggests that the kurtosis n remains finite at least up to 22 = 0.85, in agreement 
with our theoretical analysis; the numerical value of k = 3J-(0) is found to be ~ 12 for 22 = 0.85. From this table, one can 
extrapolate /1 to be r; 3.6 and A 2 (In) to be « 0.015 for 22 = 0.9. 

the distribution of returns will also display the same power-law tail. The values of fi that we find using the above fit 
are summarized in Table I. From Fig. 3, however, we see that the apparent slope of lnP(u) vs. u in the available 
range of 'large' u values is slightly smaller than the value of \i obtained from a global fit with Eq. l|3l)[l . For example, 
for z-i = 0.85, we find fx w 4, but the apparent slope is w 3.5, interestingly closer to the value reported for stocks /j«3 
Q. A slightly larger value of 22 = 0.9 would be in even better agreement with this empirical value of the exponent 
(see Table I). 

We have also tried to fit P(u) assuming an inverse Gamma distribution for the pdf of a which corresponds 

to a Student distribution for the returns. In terms of u = lncr, this reads: 

P{u) = Z' exp [-Ae- Bu - fin] , (31) 

where A and B are parameters, and fi is the power-law tail of the return distribution. Of course, this distribution 
cannot be exact in the present case since it takes non zero values when u < 0. Although it is definitely not as good a 
fit as Eq. (|30|l . it is quite acceptable, meaning that returns are indeed close to being Student distributed in our model. 
On the other hand, a log-normal distribution for a is clearly inadequate to describe our data (it would correspond to 
a parabola in Fig. 3.) 

From Table I, we see that (a) the numerical value of the average volatility is close to its theoretical value up to 
22 ~ 0.75, beyond which a systematic underestimation of the true volatility a 1 is observed, which reaches 10% for 
22 = 0.85; (b) the kurtosis n — 3^(0) increases with 22, as expected, and seems to remain finite at least up to 
22 = 0.85, beyond which /i appears to drop below 4, signaling a divergence of k (and correspondingly an even more 
difficult determination of the statistical properties of the system). 
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FIG. 4: Left: variogram of the square volatility for 22 = 0.85, a = 1.15 and a = 1.3 and fit with power-laws with v = 2(a — 1). 
We also show the data for US stocks with one day corresponding to 300r for a — 1.15 and 80r for a — 1.3, for which the 
agreement is clearly not as good. Right: variogram of the log-volatility for 2:2 = 0.85, a = 1.15, and fit with an affine function 
of lni?, the slope of which yielding (twice) the intermittency parameter A 2 , here found to be » 0.0125, whereas the US data 
suggests a larger value A 2 ~ 0.0165 - which would be matched by choosing 22 = 0.90, for which we estimate from Table I 
A 2 w 0.015. 



2. Volatility correlations 



We now turn to the temporal correlations of the volatility. Several characterizations of the "long-memory" property 
are interesting to consider. Well studied quantities are correlations of different powers of the volatility, or of the 
logarithm of the volatility. In our model, we of course know exactly the volatility at any instant of time, whereas, as 
pointed out above, in real conditions one only has access to price changes, from which a (noisy) proxy of the volatility 
is constructed. We find numerically that the shape of the correlation function can be noticeably different for these 
two quantities when the noise is large; this observation may be especially important for calibration. 

From Table I, one sees that the average volatility is rather ill-determined in the cases most relevant for applications, 
i.e. Z2 and a both close to unity, variograms should be preferred to correlograms |ll). In other words, we will study 
the following quantity: 

V n {t) = - a^tf). (32) 

These variograms are plotted in Figs. 4 a,b for the case n — 2 and n — ► 0. This last case reproduces, thanks to 
the 1 /n 2 normalization, the variogram of the logarithm of the volatility which has been much studied in the context 
of multifractal models (see below). The case n = 2 is important because it can be analytically studied, as we did 
in Section II I II and because it is related to the kurtosis of the distribution of returns for different time lags, which 
determines the smile of option prices. Form Eq. 1|26[) . one finds that for large t, one should observe: 

V 2 (£) ~ 2^o - 2Foof - 2^r 2l/ + ... (33) 

When a is close to 1, v is small and one should a priori be prepared to see corrections to the asymptotic result coming 
from the l~ 2v contribution. Fig. 4-a however shows that, for 22 = 0.85 and v = 0.3 our numerical result is rather well 
fitted by the dominant term of Eq. I|33|l . The value of the apparent exponent v however increases when zi decreases 
(in which case the contribution of the subleading term becomes more important). We also show the US stock data 
(that corresponds to v ss 0.25 ^l| ), which can be matched quite well with the model. In order to test the sensitivity 
of V2[l) to the value of a, we also show in Fig. 4-a the case a = 1.3, corresponding to v = 0.6. The agreement with 
empirical data is clearly not as good, a conclusion confirmed by all other observables we studied. 

Another interesting quantity, less noisy than the square volatility, corresponds to n —t 0. As discussed below, this 
log-variogram appears naturally in the context of multifractal models. The result for n = is shown in Fig. 4-b; we 
see that it can be fitted approximately by the multifractal prediction [30j : 



VoW~2A 2 ln|, 



(34) 
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Gaussian limit 

z,=0.7, 0.75, 




z,=0.70, 0.75, 0.80, 0.85 
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FIG. 5: Evolution of two cumulants of the distribution of returns with the lag I, for different values of Z2. Left: rescaled mean 
absolute deviation T(£). Note that the evolution is non- monotonous as a function of I. Right: Excess kurtosis k(£). We have 
also shown (symbols) the corresponding cumulants for US stocks, where we choose I = 300 to correspond to one trading day, 
as in Fig. 4. 

where A 2 is called the intermittency parameter, and T = Lt is usually called the integral time. The value of A 2 for 
different values of z-i is given in Table I, together with another determination of A 2 discussed below. For Z2 — 0.85, 
a = 1.15, we find A 2 « 0.0125, whereas our US data gives A 2 » 0.0165, or A 2 » 0.018 for the S&P100 Index, given 
in [53. This suggests that the optimal value of z-i might in fact be closer to 0.9, for which we estimate from Table I 
A 2 « 0.015. This conclusion is reinforced by the analysis of section VI. 

B. Distribution of returns over different time scales 

Since the noise variable £ is Gaussian, one can obtain the distribution of returns on the elementary time scale £ = 1 
from the distribution of the instantaneous volatility a. For example, an inverse Gamma distribution for a leads to a 
Student-Tsallis distribution for r. As discussed above, the actual distribution of volatility in our model appears to be 
slightly different from an inverse Gamma distribution; therefore the distribution of returns in our model will be close 
to, but different from, a Student distribution. On larger time scales, the distribution progressively becomes Gaussian. 
However, the convergence is very slow precisely because of the long-memory of the volatility, parameterized by the 
exponent v. A way to quantify this convergence is to measure the cumulants of the distribution, for example the 
excess kurtosis k(£), expected from our theoretical analysis to decay as t~ v , or the rescaled mean absolute deviation 
T (£), defined as: 



T{£) = 



1 



\ x i+t ~ x i 



(35) 



For a Gaussian distribution, one should find T = These quantities are plotted as a function of £ in Figs. 5-a,b, 

for different values of z-i and for v = 0.3. An a priori unexpected feature is that non-Gaussian effects actually first 
increase for small £, before decaying back to zero beyond a certain £ — £* 50. The origin of this non monotonicity 
was discussed in section UTTl and is clearly related to the assumption that the noise £i is Gaussian. Any extra kurtosis 
coming from unpredictable jumps in the price, not captured by the feedback mechanism of our model, will strongly 
affect the shape of T(£) and k(£) on short time scales, and remove this non-monotonicity which, to the best of our 
knowledge, is not observed on empirical data, even on very short time scales. Another possibility is to change the 
shape of ge for small ts. 

Of course, the knowledge of n and T is insufficient to fully characterize the whole distribution on different time 
scales. We have in fact found that a Student-Tsallis distribution with a time dependent number of degrees of freedom 
is an acceptable fit of this distribution for all values of £. In line the notation of ref. |4l| 1 we write this distribution 



A 



(3- 5 )/(5-l) 




(A 2 + (,7-l)A 2 )i/(9-i) 



(36) 
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with an £ dependent parameter q(£). In the limit q — > 1, the distribution becomes Gaussian. If the distribution is 
indeed given by Eq. I|36[) . the relation between T and q reads: 



/ 2(^-2) r(^) _ 6 

with /x = (3 — q)/(q — 1). Using these relations, one can infer, from Fig. 5-a, the value of q that one should use for 
different times scales in order to get an approximate functional form for the distribution of returns. This is useful for 
option pricing, for example |4lj . 



C. Multifractality 

A property related to the systematic change of the distribution of returns with I is multifractality, which means 
that different moments of price changes scale as a power of time, but with different scaling exponents. More precisely, 
multifractal scaling is the following property: 

M n (e) = (\x i+i -Xi\ n )=A n £^; £^L (38) 

where A n are constants and £„ is an n-dependent exponent. In the monofractal case, where the distribution is the same 
on all time scales up to a rescaling of the returns, then = njlC^. The simplest example is obviously the (geometric) 
Brownian motion, for which Q n = n/2. Any deviation from a linear behaviour of £ n is coined multifractality, for which 
several explicit models were proposed recently EI 

One example is the Bacry-Muzy-Delour (bmd) stochastic volatility model, which makes the following assumptions 

MM 

• the log- volatilities lner^ are multivariate Gaussian variables (or more generally infinitely divisible |3l|). 

• the log-volatility variogram is given by Eq. 1341) 

• the volatilities o~i are independent from the (Gaussian) noises 

From these assumptions, one can compute exactly the moments of the return distribution on different time scales. 
One finds that these are indeed given by Eq. (|38(l . with: 

Cn = ^[l~A 2 (n-2)], (39) 

whenever n < 1/A 2 , beyond which the moments are infinite. (All A n 's can also be exactly computed These 
assumptions and predictions were found to account rather well for some aspects of empirical data. 

In the present section, we show that although our model is, strictly speaking, not multifractal, many of the mul- 
tifractal predictions actually hold numerically quite accurately. This means that our model can account very well 
for apparent multifractal properties of financial time series, and in fact cures some of the deficiencies of standard 
multifractal models (see below). First, we note that our model is not multifractal since the moments M n (£) can 
be exactly computed to be sums of power-laws with different exponents, and not a unique power-law (see 25] for a 
related discussion). For example, M\{€) is the sum of £ 2 , £ 2 ~ v , ( 2 ~ 2v ; etc., and therefore does not scale as a unique 
power-law. However, as we show now, the numerical behaviour appears difficult to distinguish from a unique, effective 
power-law. j^] 

We have computed numerically M n (£) for £ > £*, where £* corresponds to the maximum of k(£) or the minimum of 
T(£) appearing in Figs. 5-a,b. In this regime, one can neglect the contribution of terms involving (er 2 £ 2 ) c , and M n {£) 
can be expressed as: 

(t \ "/ 2 
PiV h (40) 

which is the quantity that we studied numerically, because it is much less noisy than the direct calculation of moments 
of returns. The results are shown in Fig. 6-a in a log-log representation, for z 2 — 0.85 and a — 1.15, from which it 
is obvious that pure power-laws arc indeed excellent fits. From the slope of these lines one obtains the exponents £„, 
shown for different values of z-i in Fig 6-b. We note that: 
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FIG. 6: Left: Evolution of different moments M n (£) of our model with £ in a log-log representation, which allows one to extract 
from the slope of these lines, the exponents shown in the inset. Also shown in the inset is the parabolic fit suggested by the 
BMD model, Eq. (1391 . Right: The exponents £ n as a function of n for different values of the feedback parameter Z2, and the 
corresponding results of US stocks (triangles) . 



• A parabolic fit of £ n as a function of n, as in Eq. (|39fl . gives an excellent representation of our data (see Fig 
6-a, inset). 

• From this parabolic fit, a value of A 2 can be extracted for different values of z^- This value of A 2 is four times 
larger than that extracted from the variogram of the log- volatility, in contradiction with the BMD model, where 
both should be equal (see Table I) . However, we note that a similar discrepancy with the BMD model is observed 
on US stock data as well, but that both observables are fully compatible with our model with the same set of 
parameters. The discrepancy with the log-normal BMD model is due to the underestimation of the probability 
of large events in that model. 

• The intermittency parameter A 2 , that gauges the degree of multifractality (i.e. the deviation of £ n from a 
straight line), increases as Z2 increases. The multifractal spectrum extracted from US stock data corresponds to 
A 2 ~ 0.055 and matches quite well our numerical points for Z2 — 0.85. Similar values of A 2 have been reported 
for other markets as well (see e.g. [24l l30| ). 

The bmd multifractal model makes other, even more detailed predictions, about the relaxation of volatility after a 
volatility shock. We now turn to this topic to show that our model can also reproduce these more subtle features. 



D. Response to volatility shocks 

A question of great importance for option pricing and risk management concerns 'aftershocks'. It is well known 
that after a large market move, the volatility remains high for a while. The precise question therefore is: conditioned 
to a large volatility burst, how fast will the market revert back to normal? This has been addressed both empirically 
and theoretically, within the context of the BMD model [53J. One finds that after a shock, the volatility reverts to its 
normal level very slowly, as a power-law of the time i after the shock: 

Aa l+e ~ A £- e , (41) 

where i is the time of the initial shock, Ac the excess volatility over its average value, and Ao the amplitude of the 
initial shock. For rather large shocks, the empirical data suggests w 1/2 |53l l5^ |. while decreases for smaller 
shocks. Interestingly, the multifractal BMD model suggests that the exponent 9 in fact depends continuously on the 
amplitude of the initial shock, and decreases from the value 1/2 as the amplitude of the shock decreases j53j. This 
prediction was found to be in remarkable agreement with empirical findings, giving strong support to the bmd picture. 

We have therefore computed the volatility relaxation process within our model, following the methodology of |53| . 
We compute the average volatility a time I after the shock, conditioned to a shock of a certain amplitude. The 
relaxation curves are shown in Fig 7, again in the case z-i = 0.85, a = 1.15. We observe that the predictions of the 
BMD model are again quite accurately verified by our model, which, by the same token, is an alternative candidate to 
explain empirical results. 
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FIG. 7: Evolution of <r 2 conditioned to an initial volatility (<r 2 )e 2s , with s 
lines are power-law fits, with exponents 0(s = —1/2) « 0.22, 0(s = 1/2) ~ 
quoted in for the S&P100. 



= -1/2, 0, 1/2, 1 (from bottom to top). The dashed 
0.17 and 9(s = 1) « 0.30, very similar to the results 



Another characterization of aftershocks inspired from research on earthquakes is the so-called Omori law, which 
states that the probability of an aftershock larger than a certain threshold occurring a time i after the main shock 
decays as l/£ p , with p as 1. This law was checked for stock markets in [5J| on a handful of 'significant' crashes. In 
our model, crashes are self-generated and not related to external news, obviously absent from the model. We show 
in Fig. 8 the numerically determined Omori law after large endogenous crashes, for which we obtain a significantly 
different value of p w 0.5, compatible with the value of 8 reported above for large crashes. On the other hand, it is 
easy to compute the volatility response to an exogenous crash, represented by a large instantaneous jump added 'by 
hand' in the time series. If the amplitude of the jump at time i = is J, the volatility after the crash is given by: 

((t|)jw (a 2 )+gJ 2 M e . (42) 

Using Mi ~ £~ a , we find that the probability of an aftershock larger than a given threshold also decays as £~ a for 
large enough t. Since the value of a is close to unity, an approximate Omori law with p w 1 will be observed after 
anomalously large, exogenous crashes in our model. The distinction between endogenous and exogenous crashes, 
suggested in |53|. makes perfect sense in the context of the present model, where endogenous crashes are, in a precise 
sense, the result of progressive volatility built up, resulting from the ARCH like feedback effect. This volatility built 
up is in fact related to the non monotonous behaviour of the kurtosis in our model. 



E. Time reversal symmetry 



A question of general interest is whether financial time series 'know' about the arrow of time, i.e. whether it 
is possible to compute any observable that distinguishes past from future (see |5(|). Although the answer of this 
question would appear, to the layman, to be trivially yes, things turn out to be much more subtle, and of considerable 
importance. For example, the usual Brownian motion, all Levy processes and all multifractal models constructed up 
to now (including Mandelbrot's cascade, the BMD model or the version studied by Lux in [^H) are strictly invariant 
under time reversal symmetry (trs)! Financial data, on the other hand, do reveal non trs effects. A simple example, 
on which we will expand in the next section is the leverage effect, which is a causal correlation between past price 
changes and future volatilities: a drop in price leads to an increased volatility. This effect in turn leads to some 
(negative) skewness in the distribution of returns (see below). 

Here, we want to discuss a distinct effect, recently evidenced by Zumbach and Lynch In order not to mix 

this effect with leverage, one can study FX rates between two large currencies, for example Euro vs. Dollar. In this 
case, any leverage correlation or skewness, if present, is very small. In spite of this, there is a clear time asymmetry 
in the volatility process: as shown in |13| . the correlation between large scale, past volatilities and small scale future 
volatilities is larger than between small scale, past volatilities and large scale future volatilities. This effect was 
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FIG. 8: This Omori plot shows the cumulative number of aftershocks (i.e. returns with an amplitude larger than a certain 
threshold) following a main shock, and fit with N(£) ~ Main shocks were defined as returns larger than a third of the 
maximum return observed over the whole time series (of length 850,000), and aftershocks as returns larger than a third of the 
main shock. 
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FIG. 9: Zumbach's mugshot for our model: contour plot of the correlation between past volatility and future volatility, measured 
on different time scales. For a trs process, this mug-shot would appear symmetric around the diagonal, whereas empirical 
data shows, as in this figure, that the region below the diagonal carries more correlation than the region above it. 

also noted in 57], but on the example of the S&P 500 index for which the leverage effect is very strong. We have 
computed this correlation in our model, following the methodology outlined in [l^, where the idea of 'mug-shots' was 
introduced to represent graphically such past volatilities/future volatilities correlations. The mug-shot corresponding 
to our model is shown in Fig. 9. It is clear that our model - almost by construction - captures such a non trs effect. 
This was already noted in [13( for a similar model. 

We think that the time asymmetry revealed by Zumbach's mug-shots is extremely important: first, it imposes a 
theoretical constraint on the eligible models of financial time series that most of them fail to obey. Second, it is a 
direct proof of the existence of feedback effects in financial markets: the history of past price changes does have a 
direct impact on the decision and behaviour of traders - in plain contradiction with the efficient markets dogma. 
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V. THE LEVERAGE EFFECT 



Up to now, we have only discussed our feedback model under the assumption of a symmetric reaction of the market 
participants to price changes. If negative price changes have a larger impact than positive price changes, i.e, if the 
distribution of thresholds shown in Fig. 1 has some asymmetry, one will observe negative correlations between past 
price changes and future volatilities (leverage effect) and some skewness in the distribution of returns, totally absent 
from the above model. The natural way to generalize Eq. @ to account for such an asymmetry is to write: 



l) 



1=1 



(Xi - Xi-i) 2 



(43) 



where ip measures the strength of the asymmetry. The case ip = reproduces the model studied above, while ip < 
induces a leverage effect. A sufficient condition on ip that ensures that a 2 always remains positive is to impose that 
each term of the sum over I contributes positively. Writing as an identity 1 = X)^=i5^/ Z 2, one obtains: 



<pX H > 

Z2 



VA, 



(44) 



or: z 2 (p 2 < 4. 

The leverage correlation can be defined as:j 



r 2\3/2 



(45) 



This quantity is found empirically to be close to zero for i < j and negative for i > j. It is not difficult to compute 
exactly this correlation function in our model, provided the distribution of ^ is even. 
One finds: 



(46) 



One should also note that the average volatility is unchanged by the leverage term ip. Therefore, using (a 2 ) = 
tTg/(l — z%), we finally find: 



00 1 

C(£) = pgVT^Y, -TH^ £l/2 ~ a - 



(47) 



The decay of the empirical leverage correlation with lag, although noisy, can be fitted by a power-law of exponent 
close to 0.5, not far from a — 1/2 (see Fig. 10). A power-law decay of the leverage correlation was also proposed in 
the context of the multifractal BMD model in [33, y4| . 

This quantity is important because it governs the behaviour of the skewness of the return distribution on different 
time scales. It is indeed easy to show that the normalized skewness of returns on scale £, S(£) is given by [Tl| : 



(48) 



From the above expression, we see that even if the return distribution is symmetric on the smallest time scale 
(5(1) = 0), a negative skewness appears for I > 1 when ip < 0, and decays back to zero for very large lags. However, 
once again, the proximity of the critical line a = 1 beyond which the process becomes non stationary, leads to a very 
slow decay of the skewness. Empirically, for daily returns of individual stocks, one finds S w —0.1, corresponding to 
tp ~ —1 when a = 1.15 and Z2 — 0.85. 

The skewness of stock indices, on the other hand, is generally much larger (by a factor 10) than that of individual 
stocks. This is due to an enhanced downside correlation, which should be modeled using the multi-asset model 
discussed below. 

Note that the extra asymmetric term introduced in this section actually contributes also to the volatility-volatility 
correlation T computed above and also to the kurtosis. For large £, this extra contribution behaves as (1 — Z2) l p 2 /f 
and adds to the dominant term computed in section ITTT1 
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FIG. 10: The leverage correlation C (with a minus sign) as a function of lag, for US stocks, in a log-log representation. The 
straight line corresponds to the best power-law fit over the whole range and has slope —1/2, whereas the prediction of our 
model for a = 1.15 is a slope of —0.65 (dashed line). Note however the scatter in the empirical data points: the leverage effect 
for stocks is weak and hard to measure [8|. Inset: same data in a linear representation, with the prediction of our model. 



VI. SOFT CALIBRATION WITH REAL TIME SERIES 



We have shown in the above sections, using both analytical arguments and numerical simulations, that our model 
Eq. © is able to reproduce semi-quantitatively many of the stylized facts of financial time series that have been 
reported and studied in the literature. We have in fact shown, in many of the above figures, empirical data that 
match quite well, at least to the eye, the predictions of the model. What do we mean by 'semi-quantitatively'? Can 
one be more quantitative and calibrate, in a standard econometric sense, our model to empirical data? 

We believe that our model is interesting precisely because it clearly underlines the limits of such an ambition. The 
empirical data clearly suggests that any faithful statistical model of financial time series must be somehow close to 
being non-stationary. This is obvious from the very existence of option markets, which demonstrate the difficulty of 
measuring and predicting the volatility, even on rather long periods: the at-the-money vol of long-dated options still 
moves around quite a bit from day to day and there is a persistent smile, symptomatic of a long-memory extending 
to a few years |l ll . We have shown in the above figures that empirical data on stocks seem to favor values of Z2 and a 
that drive our model very close to instability. This means that even a million step long simulation of our theoretical 
model for realistic parameters is insufficient to determine the true value of the volatility to better than 5% (see Table 
I); such an uncertainty affects all the observables of the model. How can one believe that anything more precise than 
this can be reached on real empirical data? Available time spans are necessarily restricted, true jumps and overnight 
effects make the returns even more kurtic, true seasonalities (day, week, month, quarters, years) certainly play a 
role, and non-stationarities (for example, the acceleration of the trading frequency with time) plague any attempt to 
represent the dynamics of financial markets with fixed values of the parameters on very long time scales. No test, and 
no model, should aim at more precision than reality itself. 

In this situation, we think that the only reasonable strategy is what one could call 'soft calibration', in the following 
sense: instead of focusing on a few observables that one tries to reproduce as accurately as possible to calibrate the 
model (which is always possible), one should instead find a set of parameters that approximately accounts for as 
many different observations as possible, and cross check the overall consistency of the model. This consistency is 
more important and more stringent than a perfect fit of an intrinsically elusive target. Calibration in these extreme 
conditions is an ill-posed problem that, we believe, must be supplemented by intuition on what is important and 
plausibility. 
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A. Calibration on stylized facts 



How does this work in practice for the model we studied? In the simplest version that we developed, the model 
has three important parameters: r, z 2 , a (the value of <Jq merely sets the scale of the returns, but has no bearing 
on the structural properties of the model). Accounting for the leverage effect adds one more parameter, <p. But we 
already know, both from the numerical results that show that our model leads to a non-monotonous kurtosis, and 
from common sense, that unpredictable jumps must be present and should be factored in through a non-Gaussian 
noise term £, which most probably has itself fat, power-law tails [HjllHsj. This adds at least another parameter, which 
would play an important role in an extended formulation of the model. Neglecting for now this extra complication, 
our strategy is based on the idea that different observables probe differently the influence of all parameters. This is 
if fact how we organized the numerical results of section llVl 

• The distribution of the volatility or of the returns probes primarily the value of z^. The tail exponent /i, and 
any measure of non-Gaussianity helps restricting the range of acceptable values of z 2 (see Fig. 3 and Table I). 

• The temporal correlations of the volatility is primarily sensitive to the value of a, and can be used to limit the 
acceptable range of this parameter (see Fig. 4-a), whereas the correlation of the log- volatility is most sensitive 
to the value of z 2 (see Fig. 4-b and Table 1). 

• The evolution of the non-Gaussian cumulants k and T is sensitive to Zi , a, but also to the value of the elementary 
time scale r (see Figs. 5-a,b). This has enabled us to fix r « 1/300 day to be consistent with empirical data. 
Of course this leaves us with the task of curing the unfriendly looking short scale kurtosis, but as mentioned 
above, this could be easily be dealt with a non-Gaussian noise £ . This however means that the optimal value of 
Z2 would be slightly reduced, since part of the kurtosis would already be accounted for. 

• The multifractal analysis provides a stringent cross-check of the choice of parameters, since the multifractal 
spectrum is quite sensitive to the value of z 2 (see Figs 6-a,b). 

• The consistency of the model can be probed further by analyzing the response of the volatility to shocks of 
different amplitudes, and studying the Omori plots (see Figs 7, 8). An acceptable description of this rather 
subtle statistics is, we believe, another useful constraint on the parameter range. 

• Interestingly, the leverage correlation is totally decoupled from other observables and can be determined inde- 
pendently from the study of the asymmetry of the distribution of returns and asymmetric volatility correlations, 
that allow one to fix the parameter ip. [Note however that ip ^ adds a contribution to the kurtosis of the 
process]. 

Following these steps is how we 'calibrated' our model on the average behaviour of 252 liquid US stocks in the 
four-year period 2000-2003. From Figs. 4-6, we see that the value a = 1.15 allows one to capture correct time 
dependence of the volatility correlation and of the evolution of non-Gaussian cumulants, whereas the choice of z 2 in 
the range 0.85 — 0.90 allows one to capture the correct level of non-Gaussianity and multifractality (the parameter A 2 
appearing in Table I and in Figs 6). These values of z 2 and a allow us to reproduce quite satisfactorily the whole set of 
observables that we have studied, in particular the Student-like shape of the distribution of returns with a power-law 
tail index in the right range, and the slow decay of the volatility correlation and of the kurtosis. The choice of the 
time scale r is dictated by Figs. 4-5, and is found to be on the order of l/300th of a day (a few minutes). The value 
of both Z2 and r will probably be affected by the inclusion of a non-Gaussian noise £ - we leave the detailed study of 
this effect for future work. 



There is another, more direct way to test the consistency of our model, which to some extent avoids the problem 
of the non-Gaussian nature of the noise £ (but is still confronted with the intrinsic problems of long memory and non 
stationarity) . The idea is to fix the exponent a and to regress, on empirical data, an estimate of the daily square 
volatility on the computed feedback strength, defined as: 



In practice, we have estimated a noisy proxy of square volatility of a stock as of = (H — L + \0 — C|) 2 /40 2 , where 
O, H, L, C stand for Open, High, Low, Close. We have computed Xi using the open prices and truncated the sum 



B. Volatility prediction 




(49) 
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FIG. 11: Scatter plot of of vs. Xi computed daily for 252 US stocks during the four-year period 2000-2003. The coordinates 
of each point were rescaled by the average square volatility of the stock during that time period. A moving average over 1500 
points was performed, unveiling the nearly linear average behaviour of of on Xi, assumed in our model. One can even notice 
a negative curvature for large X, as suggested by the saturation mechanism we invoked in section[n]to motivate the model. 

over £ beyond 500 days, which of course is not very accurate because when a is close to one, the above sum converges 
only very slowly. 

We then plot for all stocks of / (a 2 ) vs. Xi/(a 2 ). Using our data set, this gives ~ 400,000 points; the correlation 
coefficient between the two sets is found to be 0.285. This value is rather high in view of the roughness of our volatility 
proxy. The result is shown in Fig. 11, where we have performed a moving average over 1500 points. As one can see 
the assumption of a linear relation between of and X is rather remarkably borne out, over a rather large range of Xi. 
From the slope and intercept of the linear relation, we obtain a direct estimate of Z2, which we find to be ~ 0.9, quite 
close indeed to our previous determination. This direct estimate shows that (a) the basic assumption of the model, 
that past price changes feedback in the volatility as in Eq. (JSJ, seems to be realistic; and (b) the model is indeed 
rather close to an instability, with a feedback mechanism that leads to a substantial increase of the volatility. 

The direct determination of a using this method is however difficult: one could think of varying a and choosing the 
value corresponding to the maximal correlation between Xi and Oi . Unfortunately, the dependence of this correlation 
coefficient on a is weak and does not allow to extract a meaningful minimum, although one can see that a = 1.15 
is indeed in the range where the correlation is largest. One could also extend the above method to account for the 
leverage effect and estimate directly the asymmetry parameter ip. 

C. Summary 

In summary, we have shown that using a variety of different observables, the range of acceptable values for the 
parameters of the model can be approximately determined. We have found that using these parameters, all stylized 
facts can be quantitatively accounted for. Due to the proximity of the unstable regime, however, a very precise 
determination of optimal parameters seems illusory. On the other hand, the basic assumption of our model, that past 
price changes feedback in the volatility through Eq. @ , is rather convincingly supported by the results shown in Fig. 
11. 



VII. GENERALIZATION TO MULTIASSET MODELS 

An interesting generalization of the above model concerns the multiasset situation, for example baskets of different 
stocks with cross-correlations both in the returns and in the volatility. An obvious generalization of our model to this 
case reads: 



x i+l x i ~ r i — a i£,i ) 



(50) 
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where i denotes the time index and a labels the stocks. The £"'s are characterized by certain correlation matrix 
C a b = encoding the usual sectorial correlations. For the erf's, we write, in full generality: 

°i = <*\ 

We leave the investigation of this rich model for future work; thanks to the matrix structure of the feedback effect H 
and G, one can reproduce a large variety of volatility cross-correlations and leverage effects. Here, we note that the 
average volatilities obey the following matrix equation: 

£ U* {<7 b2 )=<, (52) 
b \ e=i a o ) 

leading to a criterion for the stability of the model, which is that the smallest eigenvalue of the matrix on the left 
hand side of this equation must remain positive, generalizing the above criterion 1 — z% < 1. 

From Ea. (|51[) one can also estimate the leverage effect for index returns, which can be much enhanced if the matrix 
H ab has large off diagonal values compared to G ab , meaning that a downward move on any other stock b is perceived 
as a source of risk for stock a, and triggers extra trades on all other stocks as well. 
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VIII. CONCLUSION AND PERSPECTIVES 



In this work, we have proposed and studied, both analytically and numerically, a multiscale feedback model of 
volatility. This ARCH-like model (similar to the one studied by Zumbach in assumes that the volatility is 

governed by the observed past price changes on different time scales, which, we argue, directly influence the activity 
of traders. Assuming a power-law distribution of the time horizon of different traders, we obtain a model that captures 
most stylized facts of financial time series: Student-like distribution of returns with a power-law tail, long-memory of 
the volatility, slow convergence of the distribution of returns towards the Gaussian distribution, multifractality and 
anomalous volatility relaxation after shocks. The model, at variance with recent multifractal models that are strictly 
time reversal invariant, reproduces the time asymmetry of financial time series revealed by Zumbach's mug-shots: 
past large scale volatility influence future small scale volatility. 

The most important conclusion of our work is the following: in order to quantitatively reproduce empirical obser- 
vations, the parameters must be chosen such that our model is 'doubly' close to an instability, i.e. two parameters are 
close to values beyond which the process becomes non stationary. This means that (a) the feedback effect is important 
and substantially increases the volatility, and (b) that the model is intrinsically difficult to calibrate because of the 
very long range nature of the correlations and the slow convergence of all observables. However, by imposing the 
consistency of the model predictions with a large set of different empirical observations, a reasonable range of the 
parameters value can be determined. Furthermore, the adequacy of the basic assumption of our model, i.e. that the 
instantaneous volatility is directly related to a power-law superposition of past square returns on different time scales, 
can be directly assessed. The model can easily be generalized to account for jumps (a feature needed to correct an 
unrealistic non monotonous behaviour of the kurtosis), skewness and multiasset correlations. 

The interest of this type of models, compared to (multifractal) stochastic volatility models, is that their fundamental 
justification, in terms of agent based strategy, is relatively direct and plausible. We believe this is a strong constraint 
which should guide the construction of any mathematical model of reality. On the other hand, our fundamental 
assumption, Eq. lO, is in contradiction with the efficient market hypothesis, which asserts that the price past history 
should have no bearing whatsoever on the behaviour of investors. The large correlation that we find between past 
price changes and present volatility (see Fig. 11) indicates that this influence is in fact quite strong. This result is, in 
our view, yet another direct piece of evidence against the efficient market hypothesis, and a clear mechanism leading 
to excess volatility in financial markets. 

Turning to financial engineering applications, such as risk control and option pricing, our model provides a well 
defined procedure to filter the series past price cha nges, an d to compute the probabilities of the different future paths. 
Similar models have been shown to fare rather well |46ll48| . Once 'softly' calibrated, the model can in principle be used 
for VaR estimates and option pricing. However, its mathematical complexity does not allow, in general, for explicit 
analytical solutions and probably one has to resort either to approximate treatments or to numerical, Monte-Carlo 
methods. The difficulty of long-memory models is that the option price must be computed conditional to the whole 
past history, which considerably complexifies both analytical solutions and Monte-Carlo methods. In other words, 
both the option price and the optimal hedge are no longer simple functions of the current price, but Junctionals of the 
whole price history. Finding operational ways to account for this history dependence seems to us a major challenge, 
on which we hope to work in the near future. 
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